Multiscale reconstruction of time series 
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A new method is proposed which allows a reconstruction of time series based on higher order 
multiscale statistics given by a hierarchical process. This method is able to model the time series 
not only on a specific scale but for a range of scales. It is possible to generate complete new time 
series, or to model the next steps for a given sequence of data. The method itself is based on the 
joint probability density which can be extracted directly from given data, thus no estimation of 
parameters is necessary. The results of this approach are shown for a real world dataset, namely for 
turbulence. The unconditional and conditional probability densities of the original and reconstructed 
time series are compared and the ability to reproduce both is demonstrated. Therefore in the case 
of Markov properties the method proposed here is able to generate artificial time series with correct 
n-point statistics. 

PACS numbers: 02.50.Ga, 05.45.Tp 
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I. INTRODUCTION 

A typical feature of complex systems, like turbulence 
and finance, is that these have hierarchical and non- 
trivial structures on different scales. Due to the hierar- 
chical structure such systems are in general very high di- 
mensional, thus there is the challenge to find appropriate 
methods to model these structures efficiently. Successful 
attempts have been made to describe these systems as 
processes in scale rather than time or space itself. Ex- 
amples are the description of roughness of surfaces JE 
turbulence [J Q earthquakes and finance 
These successful attempts are characterised by the fact, 
that they can correctly provide the joint probability den- 
sity function p(y\(r{), ■ ■■,y n {T n )) of the increments y(r) 
of the process variable x at different scales r, 



y{t,r) :=x(t)-x(t-r). 



(1) 



An increment time series y(t, r) for arbitrary scales r 
can be constructed from x(t). In the following the joint 
statistics of these increment processess are considered. 
Hereby it is assumed, that these processess, y(t,r), are 
stationary in t and the ergodic hypothesis can be applied. 
Note that the method itself is not restricted to this par- 
ticular increment definition 9| . For reasons of convenient 
notation we assume for the following discussion a given 
complex system x(t), where t denotes time but may de- 
note space as well. The smaller increments are nested 
into the larger increments. 

It should be pointed out, that the knowledge of the 
joint probability density function p(jji(ti), y n {T n )) in- 
cludes also information like multifractal scaling proper- 
ties, typical expressed by < y n (r) >~ r^ n . 

For complex systems there is a general challenge, to 
generate surrogates based on the knowledge of the sta- 
tistical properties of the system itself. One approach is 
to propose general models and proper parameter sets, 
see for example 0) . Typically remarkable features 
like correlations or the shape of single scale pdfs p{y{r)) 
are used for the parameter selection. Some of these ap- 



proaches provide also a remarkably good description of 
higher order statistics. Here we propose an alternative, 
more general approach, which is based on the knowledge 
of multiscale properties. In particular we show how a 
given time series can be extended into the future prop- 
erly. 



II. 



METHOD 



As stated above a hierarchical complex system can of- 
ten be described by means of the joint probability den- 
sity p(yi(ri), j/ n (r n )). Because of the involved scales 
the dimension of the joint probability density is very high. 
Therefore it is in general very difficult to compute it from 
empirical time series. However the description and com- 
putation can be highly simplified if Markov properties 
can be assumed. This is the case if 



p(yi{n)\yi+i(Ti-\ 



,y n {T n )) = p(yi(n)\yi+i(n+i)) (2) 



is true for all i and n > i. Without loss of generality 

we take Tj < Tj+i. It should be noted that the Markov 

property can be tested for a given data set 

this case the joint probability density can be substantially 

simplified: 



p(yi(n), -,y n (T n )) 



(3) 



p(yi(n)\y2(T 2 )) ■ ... ■ p(y n -i(T n -i)\y n (T n )) -p(y„(T n )) 

In the following we restrict the discussion to Markov 
processes in scale and to right-justified increments, i.e. 
the smaller increment is nested into the larger one and 
has the right end point in common [T^j according to Eq. 
JTJ . There are two possible ways to determine the condi- 
tional probability densities p(yj_i(Tj_i)|y,(Ti)) and the 
unconditional probability density p(y n (r n )). The first 
one is straightforward. A typical time series from the 
complex system is examined, and the probability den- 
sities are determined, by simply calculating histograms 
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from the data. Due to the finite length of the data set, 
very small probabilities tend to be estimated as zero. In 
order to obtain the joint probability, the conditional and 
unconditional probability function have to be multiplied. 
Therefore the joint probability will also be zero for this 
value. Due to that, the probability density can be under- 
estimated, especially in case of large negative and posi- 
tive increments. 

A second possibility to estimate these probabilities is 
to use the ansatz presented in 0, 0, 0, El ■ The prob- 
ability densities are now obtained as a solution of the 
Fokker-Planck equation: 



dp(y, t) 
dr 



.—D^Hv t) + —D^(v t) 
dy ' dy 2 [V ' ' 



P(v,r). (4) 



Note that this equation holds also for the conditional 
probability densities. The Kramers-Moyal coefficients 
are defined as 



(5) 



:!At 



(y' -yYp(y'(T + AT)\y(T))dy' 



and can again be directly computed from given data. The 
problem of zero values vanishes, or at least is shifted to 
increments with a very large magnitude (see for example 

0). 

Next we focus on the procedure to create time se- 
ries. We start with the most general case, that the 
time under consideration for the new element of the se- 
ries (x(t)) is named t* , then the knowledge of all points 
{x(t* — Ti), ...,x(t* — r n )} and the corresponding proba- 
bility densities for {y(ri), j/(r„)} are needed in order 
to choose the value of x(t*) correctly. Let's assume that 
x(t*) is chosen as a value for x(t*). The corresponding in- 
crements yi(ri), with their common right endpoint x(t*), 
are defined through 



Viin) := x(t*) - x(t* - T.i). 



(6) 



Because P(y n (T n )) and x(t* — r„) are known, the prob- 
ability P(x(t*)\x(t* — t„)) of this event is also known. 
But this conditional probability includes only the infor- 
mation of one scale, namely r„. To add the informa- 
tion available on the scale r„_i, p(y„_i(r„_i)|?/„(T n )) 
is also needed. These quantities then determine the 
probability P(y„_i(T n _i), y n (T„)) and contain the ad- 
ditional information of x(t* — r n _i). Therefore now 
P(x(t*)\x(t* - T n ),x(t* - r n _i)) is known. Due to the 
Markov property this can easily be extended to the prob- 
ability P(x(t*)\x(t* -T n ),...,x(t* -Tl)). 

Thus the conditional probability density conditioned 
on many different scales is known. For a finite range of 
values for x(t*) and a finite number of necessary scales 
the probability density p(x(t*)\x(t* — r„), ...,x(t* — t±)) 
is obtained. p(x(t*)\x(t* — r„), x(t* — t±)) now con- 
tains all relevant statistical information of the previous 
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FIG. 1: Three runs of reconstructed time series, which have 
the same initial conditions. Left to the vertical line, the 200 
data points, which were taken as initial condition, are shown. 
On the right side the following 200 points from our recon- 
struction are shown. x(t) has been normalised to zero mean 
and is presented in units of the standard deviation (Too of the 
whole dataset. 



time series for a correct choice of the value x(t*). Choos- 
ing now a random value from this distribution, the time 
series will be extended correctly by another point. Re- 
peating this procedure will produce a time series, which 
exhibits the correct joint probability density function for 
all considered scales t\, ...,t„. Obviously the estimation 
of p(x(t*)\x(t* — T n ),...,x(t* — ri)) has to be repeated 
for each new value of the time series. Some open ques- 
tions remain: 1. How many and which different scales 
should be included? 2. How to cope with unavoidable 
discreteness of the statistics? 3. Which initial conditions 
should be used for the simulation? The first question 
of correct choice of scales gets simpler, if for the com- 
plex system there exists a large scale, L, above which 
neighboring increments become statistical independent. 
A further simplification is given, if there is a small scale, 
the so called Markov-Einstein coherence length, l mar [l7| . 
below which the Markov properties are no longer fulfilled. 
For this case we suggest as number N of scales: 



N = 



log 2 



(J) 



In some cases, L can be very large, which is especially 
true for bounded time series. Here we propose to restrict 
to the small and medium size increments. The influ- 
ence of the large size increments can be included by a 
specific condition on the largest conditional increment. 
We close the sequence of conditional pdfs as follows: 

p{vi{n)\y2{T2))-----p{yn-i{Tn-i)\yn{Tn))-p{yn{j n )\x(t*- 

Tn)) 
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FIG. 2: a) The pdfs from the original (solid black line) and the reconstructed (circles) data set are shown. The considered 
scales are 2" • l maT with n = 7,6,5,4,3,2,1,0 from top to bottom. The pdfs for different scales are shifted for clarity of 
presentation, b+c) Contour plots of conditional probability densities p($/i-i(Ti_i)|yi(Tj)) with n = 2 • Tj_i and n-i = lmar (b) 
and Ti-i = 2 6 ■ l mar (c). The ones from original data are denoted as thin solid lines, while the ones from reconstructed data as 
thick dotted lines. Increments are given in units of the standard deviation <7oo of the whole dataset x(t). 



III. RESULTS FOR TURBULENCE 

As a real world example turbulence is chosen, where 
hierarchical Markov processes were well-demonstrated 
0, 0| . The dataset was produced by a free jet of air 
into air. The dataset consists of 12.5 • 10 6 data points, for 
more details the reader is refered to ^3 • A discretisation 
of the scale and the velocity is chosen as follows. For the 
smallest scale we choose the above mentioned Markov- 
Einstein coherence length, Z mar , in our case 6mm. The 
selected larger scales are 2™ • l mar - So a reconstruction 
of the hierarchical process is attempted for the scale in- 
terval lmar to 128 • l m ar, where the largest scale is well 
above the so called turbulent integral length. The influ- 
ence of the larger scales is included by the approximation 
mentioned above. For this choice of scales we are sure, 
that the necessary Markov properties are always fulfilled. 
Due to the discretisation of the velocity values, the ve- 
locity increments on the smallest scale l m ar can have 201 
possible values. For the reconstruction in this section, all 
pdfs were calculated directly from the data without the 
use of Eqs. (3J and (J5J. 

In order to show typical reconstructed time series, 200 
data points from the dataset described above were taken 
as an initial condition for the calculation of 200 synthetic 
data points. The results for three runs are shown in Fig. 
^ Clearly, the three different runs yield similiar results 
for the first synthetic data points, but become more and 
more different for the further data points. However the 
general characteristics of the initial time series are re- 
tained. In order to show, how well the multiscale be- 
haviour is captured, a reconstructed time series which 
contains 10 6 data points is produced, without using the 



original data set as initial condition. As an initial con- 
dition in this case, a sequence of zeros was used and 
the first elements of the simulation were discarded. The 
corresponding probability density functions for the eight 
different scales are shown in Fig. |2Ji and compared to 
those of the original data. It is clearly seen, that for 
| y | iS 2a 'oo the agreement between the original and re- 
constructed pdfs is quite good. For extreme values the 
reconstructed data deviate due to the reasons specified 
above. With increasing number of given original data the 
range of agreement in the pdfs extends. 

In order to see if the multiscale description is correct, 
not only the unconditional probability densities have to 
be examined, but according to Eq. © also the con- 
ditional probability density of first order. In Fig. 
and Fig. the results are shown for the smallest and 
the largest scales considered here. In both cases again 
the empirical conditional probabilities are modelled very 
well. 

In a next step we compare our method with a one scale 
method to reconstruct time series. Therefore a random 
walk is constructed, such that the distribution of the in- 
crements on the smallest scale is identical to that of the 
empirical time series, i.e. only p(yi(ri)) is used. The re- 
sults for this adapted random walk are shown in Fig. [3] 
As can be seen for small scales the adapted random walk 
provides good approximations for the increment distri- 
bution, but for the larger scales, the approximation with 
the adapted random walk is very poor. This illustrate 
the advantage of using multiscale information for the re- 
construction. 
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FIG. 3: Pdfs in analogy to Fig. Et are shown. For the 
pdfs (circles) of the reconstructed data now a simple one scale 
process was used. Further parameters are chosen as in Fig. 
Ill 

IV. CONCLUSIONS 

A new method has been presented, which allows for 
a multiscale reconstruction of a time series. More scales 
are modelled simultaneously, using conditional probabili- 
ties of properly chosen increments. A necessary feature of 
the underlying complex system is that it exhibits Markov 
properties of the scale dependent process. This can be 
tested statistically. In this case all conditional probabili- 



ties of arbitrary order can be substituted by conditional 
probabilities of first order, which can be calculated easily 
from empirical time series. Using the empirical probabil- 
ity densities a time series can be reconstructed. Exam- 
ples of reconstructed time series for a real world example 
from turbulence have been shown. The reconstructed 
time series reproduces the unconditional and the condi- 
tional probability densities very well. The advantage of 
this method is the use of the distribution which are ob- 
tained directly from the data. Therefore a modelling of 
an asymmetric or non-Gaussian distribution is possible. 
The reconstruction can be applied to any range of scales 
without a change of the underlying method, as long as 
Markov properties in scale are present. Furthermore us- 
ing this approach the full probability density for the next 
potential step of a time series is provided. Using an em- 
pirical time series as a starting point, this procedure can 
be applied for prediction purposes. But contrary to many 
other methods not only the mean value for the next step 
is provided, but the full distribution. This therefore al- 
lows also the prediction of volatility or the probability 
of extreme values. Due to the fact, that only empirical 
distributions are used, there is no necessity of parameter 
estimation. Also this method is quite fast, in order to 
simulate a process with 10 6 points it takes around one 
hour on a normal personal computer. In our opinion 
further improvement can be obtained, if refined uncon- 
ditional and conditional pdfs, which provide better ap- 
proximation of the real distribution for extreme values, 
are used. 

We are very grateful to David Kleinhans, Stephan 
Barth and Carola von Saldern for interesting discussions 
and remarks. 
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